#经度最大值117.0083，经度最小值115.0083
#纬度最大值11.0068，纬度最小值9.0045
import matplotlib.pyplot as plt
import numpy as np
import os
import pandas as pd
from scipy.interpolate import griddata
from matplotlib.ticker import MultipleLocator
os.chdir(r'D:\pythonProject')
plt.rcParams['font.sans-serif'] = ['Times New Roman']  # 字体
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号
plt.rcParams['font.size'] =10

df = pd.read_excel('B_x.xlsx')

s = 1/60

x = df['经度']
y = df['纬度']
z = df['重力异常值']

X = np.arange(x.min(),x.max()+s,s)
Y = np.arange(y.min(),y.max(),s)
X,Y = np.meshgrid(X,Y)
print(X)
new_Z =griddata( (x,y),z, (X,Y) ,method='cubic')# new_Z是插值后矩阵
print(new_Z)

#将矩阵导入到excal
df = pd.DataFrame(new_Z)
df.to_excel("附件一梯度异常值.xlsx", index=False)

# # 将矩阵重新形状为11x11的小矩阵块
small_matrices = new_Z.reshape(11, 11, 11, 11).swapaxes(1, 2).reshape(121, 121)#
# ##创建一个空列表来填充重力异常值的差值
my_list = []

# 输出小矩阵块
for i in range(11):
    for j in range(11):
        print(f"小矩阵 ({i+1},{j+1}):")
        print(small_matrices[i*11:(i+1)*11, j*11:(j+1)*11])
        ##计算重力异常值的差值
        min_value = np.min(small_matrices[i*11:(i+1)*11, j*11:(j+1)*11])
        max_value = np.max(small_matrices[i*11:(i+1)*11, j*11:(j+1)*11])
        chazhi_calue = max_value - min_value
        print(chazhi_calue)
        my_list.append(chazhi_calue)
        print()

# 创建一个DataFrame对象
df = pd.DataFrame({'Value': my_list})

# 将DataFrame写入Excel文件
df.to_excel('附件一重力梯度.xlsx', index=False)

fig,ax = plt.subplots(figsize=(6,5))
p = ax.contourf(X,Y,new_Z,levels=100,cmap='RdBu_r')

ax.xaxis.set_minor_locator(MultipleLocator(s))
ax.yaxis.set_minor_locator(MultipleLocator(s))
ax.xaxis.set_major_locator(MultipleLocator(0.25))
ax.yaxis.set_major_locator(MultipleLocator(0.25))

ax.set_xlabel('经度',fontproperties='Simsun')
ax.set_ylabel('纬度',fontproperties='Simsun')
ax.set_title('重力异常值插值热力图',fontproperties='Simsun')
ax.grid(which='both',lw=0.4,color='k')
ax.tick_params(direction='in',which='minor',length=0)
cbar = plt.colorbar(p)


cbar.ax.set_ylabel('重力异常值',fontproperties='Simsun')
plt.show()





